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INTRODUCTION 

With  the  advent  of  high-speed  electronic  digital  computers  it  became 
feasible  to  solve  large  systems  of  equations.  Many  techniques  have  been 
developed  to  solve  these  equations.  One  of  these  involves  calculating  the 
inverse  of  the  coefficient  matrix.  This  report  will  consider  eight  different 
methods  of  calculating  the  inverse  of  a  matrix. 

Before  proceeding  with  the  discussion  of  matrix  inversion,  it  will  be 
necessary  to  define  several  terms.  A  matrix  can  be  defined  as  a  rectangular 
array  of  elements.  For  this  report,  the  matrix  will  be  square  and  its  elements 
will  all  be  either  real  numbers  or  smaller  matrices  that  may  be  rectangular. 
Capital  letters  will  be  used  to  denote  matrices,  and  small  letters  will  be 
used  to  denote  the  elements.  The  individual  elements  of  a  matrix  A  will  be 
designated  by  a^  ,  where  the  subscript  r  indicates  the  row  in  which  the  ele- 
ment is  located,  and  the  s  the  column.  For  a  shorter  notation  the  matrix  A 
will  be  denoted  by 

A  -  (a^g). 

The  size,  or  dimension,  of  a  matrix  with  n  rows  and  m  columns  will  be  de- 
noted as  n  X  m.  When  A  is  equal  to  (aj.s),  then  the  transpose  of  A,  denoted 
A' ,  is  equal  to  (agj.) .  The  diagonal  elements  of  a  square  matrix  are  those 
whose  row  and  column  subscripts  are  the  same.  The  trace  of  a  square  matrix 
A,  denoted  T(A),  is  the  sum  of  the  diagonal  elements  of  the  matrix  A,  If  a  " 
square  matrix  has  zeros  for  all  elements  not  on  the  diagonal  then  it  is 
called  a  diagonal  matrix.  If  all  the  possible  non-zero  elements  of  a  diago- 
nal matrix  are  equal,  the  matrix  is  called  a  scalar  matrix.  The  product  of 


a  scalar  c  and  a  matrix  A  is  given  by  .  ^ 

cA  "   (caj.g). 

This  means  that  the  constant  c  is  multiplied  by  every  element  of  A,  A  sub- 
matrix  of  a  matrix  A  is  an  array  formed  by  deleting  one  or  more  rows  and/or 
columns  of  A.  Partitioning  a  matrix  is  simply  dividing  the  matrix  into 
several  submatrices.  This  dividing  can  be  accomplished  by  drawing  vertical 
lines  between  some  of  the  columns  and/or  drawing  horizontal  lines  between 
some  of  the  rows*  All  matrices  partitioned  in  this  report  will  be  divided 
into  four  submatrices  by  drawing  one  vertical  line  between  two  columns  and 
one  horizontal  line  between  the  corresponding  rows.  These  submatrices  then 
become  the  elements  of*  the  matrix  and  are  denoted  by  capital  Aj.s« 

If  two  matrices  A  and  B  are  the  same  size  and  a^^   =  b^g  for  all  r  aot 
s,  then  A  ■  B,  The  sum  of  two  matrices  is  given  by 

A  +  B  -  (a^s  +bj.s). 

This  means  that  to  add  two  matrices,  it  is  necessary  to  add  corresponding 
elements  of  each  of  the  matrices.  The  set  of  all  matrices  the  same  size 
forms  a  commutative  group.  The  identity  matrix  for  matrix  addition  is  the 
matrix  Z,  which  has  zeros  for  each  of  its  elements.  The  product  of  two 
matrices  A  and  B  of  dimension  n  x  m  and  m  x  p,  respectively,  is  given  by 


m 

A.B  -  AB  -  (21  arkbks) 
k=l 


where  AB  is  of  dimension  n  x  p.  This  means  that  the  element  of  AB  in  the 


r-th  row  arxi  the  s-th  colvmn  will  have  the  above  suitmation  for  its  value. 
Matrix  mxiltiplication  is  not  always  curomutative;  however,  it  is  associative* 
The  identity  matrix,  denoted  by  I,  for  multiplication  is  a  scalar  matrix  with 
the  value  1  for  the  diagonal  elements.  When  the  matrix  A  is  multiplied  by 
itself  A* A,  it  is  denoted  by  A  .  Multiplying  A^  times  A  would  be  A^.  In 
other  words,  powers  of  a  matrix  are  similar  to  powers  of  a  number  x  in  the 
real  number  system.  The  matrix  A  is  defined  to  be  the  identity  matrix  I, 
For  some  square  matrices  there  exists  another  square  matrix  denoted  A"  such 
that 

AA-1  »  A"^A  =  I. 

This  matrix  A~  is  called  the  inverse  of  the  given  matrix,  A  square  matrix 
that  has  an  inverse  is  called  a  nonsingular  matrix.  Even  when  the  inverse 
of  a  given  matrix  exists  it  is  not  always  easy  to  fini. 

For  every  square  matrix  there  exists  a  number  of  the  system  of  the  ele- 
ments which  can  be  associated  with  it.  In  this  paper  the  number  will  always  - 
be  real.  This  number,  denoted  |A| ,  is  called  the  determinant  of  the  matrix. 
The  determinant  of  a  2  x  2  matrix  is  defined  to  be 


lAl  - 


an  ai2 
agi  a^2 


(^11  ^2  "  ^12  ^1^' 


The  determinant  of  a  square  submatrix  of  a  matrix  A  is  called  a  minor  of  A, 
If  only  the  r-th  row  and  the  s-th  column  are  deleted  to  form  the  submatrix 
the  minor  is  denoted  as  m^.^.  The  cof actor  c   of  the  element  a^.  is  given  by 
(-1)   mj.g.  The  cof  actor  differs  from  the  minor  only  by  sign  when  either  r 


or  s  (but  not  both)  is  an  odd  integer. 

For  any  square  matrix  A,  and  for  any  choice  of  r  and  s,  the  determinant 
of  A  can  be  defined  as 


n  n 

lAl    -  ^  a^i  c^  -  21  a^s  «Js 

i»l  o"l 


where  the  c^4  and  the  c .  are  the  c of actors  of  the  elements  a_.  and  a.^  of  the 
ri         js  ^-^  J« 

matrix  A,  The  first  expression  is  said  to  be  the  expansion  by  cofactors  along 
the  r-th  row  of  A,  and  the  second  is  the  expansion  by  cofactors  along  the  s-th 
column.  By  applying  this  definition  to  the  cofactors,  it  is  possible  to  eval- 
uate the  determinant  of  any  size  of  matrix  by  calculating  the  determinant  of 
a  number  of  2  x  2  subraatrices.  The  cof actor  matrix,  denoted  by  C,  can  be  de- 
termined by  replacing  every  element  of  A  by  its  corresponding  cofactor. 

Primarily,  this  report  uses  notation  as  adopted  by  Fuller  in  his  book 
Basic  Matrix  Theor;^,  (li),  where  a  much  more  conplete  discussion  of  the  pre- 
ceding definitions  can  be  found. 

Using  matrix  notation,  a  system  of  n  linear  equations  in  n  unknowns  can 
be  written  in  the  matrix  form 

.  AX  -  G, 

From  this  form  it  is  easy  to  see  hovx  the  inverse  of  the  coefficient  matrix, 

if  it  exists,  can  be  used  in  solving  a  system  of  equations.  By  ntultiplyLng 

-1  ' 

both  sides  of  the  equation  by  A 


*-,. 


IX  -  A'-'-G 
X  -  A'-'-G. 


Thus  the  solution  for  the  system  can  be  found  by  one  matrix  multiplication. 
The  matrix  A  ;d.ll  be  used  as  an  example  to  illustrate  the  different  methods 
of  matrix  inversion  in  the  next  three  sections  where 
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INVERTING  MATRICES  USING  ELEI^SNTABY  OPERATIONS 

Before  going  into  methods  of  computing  the  inverse  of  a  matrix  it  will  be 
necessary  to  discuss  what  is  meant  by  elementary  operations.  Elementary  opera- 
tions are  operations  which  can  be  performed  on  a  system  of  equations  that  will 
yield  a  new  system  of  equations  with  the  sarae  solution  as  the  original  system. 

When  these  operations  are  used  on  a  matrix,  they  are  called  the  elemen- 
tary row  operations  for  a  matrix  and  are  defined  as  follows: 

1.  Interchange  of  the  i-th  and  j-th  row  denoted  by  R(i,  j). 

2.  Multiplication  of  the  elements  of  the  i-th  row  by  the  non- 
zero constant  c,  denoted  by  cRi. 

3.  Adding  to  the  elements  of  the  i-th  row,  k  times  the  corre- 
sponding elements  of  the  j-th  row,  denoted  by  Ri  +  kRj. 

If  the  word  row  is  replaced  by  the  word  column  and  if  the  R  is  replaced 
by  C  everyvrtiere  in  the  above  definitions,  then  one  has  the  definitions  for 
elementary  column  operations. 

Of  all  the  methods  discxissed  in  this  paper  for  calculating  the  inverse 
of  a  matrix,  those  methods  which  involve  elementary  operations  are  the  most 
easily  adapted  to  machine  computation.  The  method,  (l),  involved  in  this 


section    is  the  most  general,  for  it  uses  both  elementary  row  and  column 
operations  * 

It  will  be  shown  that  if  A,  P  and  Q  are  nonsingular  matrices  such  that 
PAQ  =  I,  then 

a"-'-  =  QP.  -  ■ 

First  premultiplying  both  sides  of  the  equation  PAQ  =•  I  by  Q 

QPAQ  =  Q. 


Since  Q  is  nonsingular,  postmultiply  both  sides  by  Q~ 


QPA  =  I. 


Therefore  by  definition  of  matrix  inverse 


a""""  -  QP. 


This  method  gives  a  technique  for  computing  the  two  matrices  P  and  Q 
using  elementary  operations.  If  the  matrix  A  is  reduced  to  I  by  a  series 
of  elementary  row  and  column  operations,  then  one  obtains  a  matrix  P  by 
applying  the  same  row  operations  to  I  and  a  matrix  Q  by  applying  the  same 
column  operations  to  I.  These  two  matrices,  P  and  Q,  are  the  ones  desired 
for  calculating  the  inverse. 

In  the  actual  computation  of  A~  by  this  method,  it  is  not  necessary  to 
keep  track  of  the  operations  that  are  made.  Instead,  by  setting  up  the  matrix 


B  = 


A  I 
I  Z 


where  Z  is  the  zero  matrijc,  it  is  possible  to  operate  on  I  at  the  same  time 
A  is  being  operated  on  by  both  row  and  column  operations.  It  is  important 
to  remember  that  all  of  the  operations  must  occur  on  A  and  that  the  matrix 
B  is  used  merely  to  keep  track  of  the  result  of  the  operation  on  I,  Once  A 
has  been  reduced  to  I 


B^> 


I  P 
Q  Z 


and  P  and  Q  have  been  determined. 

One  advantage  of  this  method  is  that  it  can  be  checked  at  any  stage.     If 
one  has  obtained 


B 


V                m 

A     I 

<-> 

\   n" 

[z   z\ 

M       Z 

then  the  work  is  correct  if 


NM 


h' 


As  an  example  of  tiiis  method,  consider  the  matrix  A  given  in  the  intro- 
duction. 
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Now  using  elementary  operations 
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which  gives 
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This  value  of  A   can  be  checked  in  the  equation  AA~  "  I. 

If  no  elementary  column  operations  are  performed  on  A,  the  matrix  Q 
will  remain  equal  to  I.  Once  A  has  been  reduced  to  I,  the  inverse  of  A  can 
still  be  calculated  by 

-1       ■ 

A   -  QP. 

However,  since  Q  is  equal  to  I, 


A  similar  argument  would  show  that  when  no  row  operations  are  performed  on 
A  the  inverse  of  A  is  given  by 


-1 


Q. 


These  are  only  two  of  the  many  variations  of  the  method  of  elementary 
operations.  One  distinct  advantage  of  both  of  these  special  cases  is  that 
it  is  not  necessary  to  form  the  product  QP, 


n 


METHODS  INVOLVING  THE  DETERMINANT  OF  A  MATRIX 

Calculating  the  Inverse  Using  the  Adjoint  Matrix 

The  adjoint  of  a  matrix  A  is  given  by  adj  A  -  C  where  C  is  the  trans- 
pose of  the  cof actor  matrix.  The  following  principal,  (U),  pp.  1G3-1H, 
makes  it  possible  to  use  the  adjoint  matrix  for  finding  the  inverse: 


A  •  adj  A  -  adj  A  •  A  -  |A|I. 


If  |A|?^  0,  then  this  eq^uation  can  be  written  as 


adl_A      adiA 
^      |A|  |A|      ^      ^' 


Therefore,  by  definition  of  the  matrix  inverse 


_1      adj  A 

^  |A| 


This  method  can  be  used  for  all  matrices  whose  determinants  are  nonzero.  It 
can  be  shown  that  the  nonvanishing  of  the  determinant  of  a  matrix  is  a  neces- 
sary and  sufficient  condition  for  the  matrix  to  be  nonsingular.  Therefore  if 
a  matrix  has  an  inverse  it  can  be  found  using  this  method. 

It  is  obvious  that  once  the  adjoint  matrix  is  obtained,  finding  the  in- 
verse id.ll  not  be  difficult.  Multiplying  the  k-th  colunm  of  A  by  the  k-th 
row  of  adj  A  will  give  the  lAl .  All  that  remains  is  to  divide  the  adj  A  by 
|A|.  The  usefulness  of  this  method,  therefore,  depends  upon  the  ease  with 
which  the  adjoint  matrix  is  obtained.  In  the  case  of  a  2  x  2  matrix  the 
adjoint  matrix  can  be  written  down  immediately.  For  given 


12 


then 


Therefore 


A 


^1  ^22 


^2  "^21 


-^ 


^11 


adj  A  =  C 


^2  "^ 
"^21   ^1 


In  other  words,  to  find  the  adjoint  of  a  2  x  2  matrix,  interchange  the 
diagonal  elements  and  change  the  signs  of  the  remaining  two  elements.  How- 
ever, in  the  case  of  a  U  x  li  matrix,  it  is  necessary  to  find  the  deterroinants 
of  sixteen  3x3  matrices,  and  for  a  5  x  5  matrix  the  determinants  of  twenty- 
five  1;  X  U  matrices  are  needed.  Because  of  the  increasing  difficulty  of  ob- 
taining the  adjoint  matrix,  this  method  should  not  be  used  for  large  matrices. 

As  an  example  of  this  method  of  finding  the  inverse  consider  the  laatrix  A 
given  in  the  introduction.  Its  cof actor  matrix  is 
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A  good  method  of  checking  the  computations  for  adj  A  is  to  form  its  product 
with  A,     For  this  example  this  gives 
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J     • 


The  adjoint  matrix  is  incorrect  if  the  product  is  not  a  scalar  matrix. 

The    !Al  is  equal  to  the  value  of  the  diagonal  elements  of  the    product 
of  A  •  adj  A.     In  this  exanple, 

lAl   =  1. 


Hence 


A      •»  adj  A. 


Inversion  of  a  Matrix  using  the  Characteristic  Equation 


The  following  method,  (5),  like  the  method  using  the  adjoint  matrix,  is 


Hi 


not  advantageous  except  with  small  matrices  or  matrices  for  which  the  charac- 
teristic eqiiation  is  known.  The  basis  for  using  the  characteristic  equation 
to  find  the  inverse  is  the  Hamilton-Cayley  theorem,  which  states  that  every 
matrix  satisfies  its  characteristic  eqiiation. 

The  characteristic  equation  of  the  matrix  A  is  the  polynomial  equation 
|rl  -A  I  »  0.  The  characteristic  equation  of  the  example  matrix  A  is 


|rl  -A| 


r-1  1    +1  -1 

2  r-U   -3  0 

0  1;  r+2  -3 

-2  -1    0  T'h 
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Evaluation  of  this  determinant  gives 

|rl  -At    «  r^  -7r^  +  l8r^  -2$r  +  1  =  0 

which  Is  the  characteristic  equation  of  the  matrix  A,      ■ 

According  to  the  Hamilton-Cayley  theorem,  "satisfying"  the  characteristic 
equation  means  that  this  same  equation  with  corresponding  powers  of  A  will 
give  the  zero  matrix  when  A  is  defined  as  I,  Therefore,  in  the  example 

A^  -  7A^  +  18A^  -  25a  +  I  =  Z 


subtracting  I  from  both  sides  of  the  equation  leaves 

A^  -  7A^  +  18A^  -  25a  »  -I. 


Factoring  out  an  A  on  the  right  side  of  the  equality  gives 
■  .       (A^  -  7A^  +  18a  -  25l)A  -  -I 


'y-'i-f! 


25 


then 


(-A-^  +  7A^  -  18a  +  251) A  =  I. 


Since  the  product  on  the  ri^t  is  eqxial  to  the  identity  matrix,  the  inverse 
of  A  is  by  defixu-tion: 


a"-'-  =  (-A^  +  71?  -  18a  +  251) . 


The  inverse  of  any  non-singular  iratrix  could  be  calculated  in  a  similar 

manner  once  the  characteristic  equation  is  known.  All  that  remains  to  com- 

-1  "^       0 

pute  A   is  to  substitute  the  values  of  A-',  A^j   A  and  I  into  the  equation  ard 

add. 


A-1 


11 

0 
-36 

-21; 


9  3 

-13  -12 

-k  13 

-26  -8  - 


39 
■78 
111 

•51 


7  -Hi 

-70  li2 

98  -35 

56  U2 


-m 

-56* 

56 

77 

-56 

Ii2 

7 

98 

-18 

18 

18 

36 

-72 

-5U 

0 

72 

36 

-36 

-18 

0 

18 
0 

-72 


25 

0 

0 

0 

25 

0 

0 

0 

25 

0 

0 

0 

0 
0 
0 

25 


,-1 


25 

13 

7 

1 

•3U 

-18 

-10 

-1 

62 

33 

18 

2 

-1; 

-2 

-1 

0 
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This  value  checks  with  results  calculated  previously.  One  disadvantage  of 
this  method  is  that  there  are  no  checks  on  the  calculations  until  they  are 
completed. 

For  some  matrices  there  exists  an  equation,  which  the  matrix  satisfies, 
vrhose  degree  is  lower  than  the  characteristic  equation.  If  this  is  the  case 
then  this  minimum  equation  may  be  used  to  f ini  the  inverse  instead  of  the 
characteristic  equation. 

Matrix  Inversion  Using  the  Trace  of  a  Matrix 

The  last  method  discussed  in  this  section  is  a  method  for  computing  the 
inverse  which  involves  the  trace  of  the  matrix.  The  fact  that  it  requires  a 
greater  number  of  operations  than  most  of  the  methods  discussed  in  this  paper 
is  a  disadvantage,  but  it  has  two  definite  advantages.  First,  it  is  not  affec- 
ted by  individ;ial  peculiarities  of  different  matrices  such  as  zero  elements  on 
the  diagonal  or  determinants  of  submatrices  being  equal  to  zero.  Second,  as 
a  by-product  this  method  will  give  the  characteristic  equation  and  the  charac- 
teristic vectors  of  the  matrix;  therefore,  if  any  or  all  of  these  are  desired 
this  method  could  be  used  to  advantage. 

Let  the  n  x  n  matrix  A  have  the  following  characteristic  equation: 

.  in  -  Al  -  r^  -  p.r^-1  -  p.r^-2  .....  p   pi  .  p 

12  -^n-l    *^n 

with  characteristic  roots  being  denoted  as  r^^,  r2,   .   .   .,  r  .     Let  the  follow- 
ing series  be  equal  to  Sj^: 

^k-^l*^2*---^^n.     ' 
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1^ 
These  are  the  sum  of  the  characteristic  roots  of  A  so 


Sk=T(A^), 


Once  the  s,  are  known,  the  coefficients  of  .the  characteristic  equation  can  be 
calculated  using  Nevrton's  recursion  formula. 


Pi  '  ^1 
2P2  -  S2  -  P^s^ 


^k  "  ^k  ■  Pl\-1  -  •  •  •  -  Pk-1^1 

^Pn  =  ^n  -  Pl^n-1 Pn-1^1 

where  the  pj^  are  the  coefficients  of  the  characteristic  equation.  This  is  a 
way  of  deriving  the  characteristic  equation  in  which  it  is  first  necessary  to 
compute  the  po>;ers  of  A. 

A  method,  (3),  has  been  developed  that  eliminates  the  necessity  of  com- 
puting the  powers  of  A  and  that  will  also  give  the  inverse  of  A,  Instead  of  "" 
the  powers,  a  set  of  n  matrices  A^,  Ap,  .  .  .  ,  A  will  be  computed  using  the 
following  formulae: 

C]_  -  T(AAq)  A^  -  AAq  -  C3_I 

Cg  -  1/2  T(AA^)  Ag  =  AA^  -  Cgl 
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<^   =  1/k  T(AAj^_i)  Aj,  =  AAi^_-L  -  Cj^I 


c^  =  1/n  T(AAn_i)  An  »  AA^.i  -  c^I 


where  A-  is  defined  to  be  I.  It  will  be  shown  that 


^k  °  Pk     k  »  1,  2,  .  .  .  ,  n 


and  that 


An  =  Z. 


If  A  is  nonsingular  and  Cj^  j^  0,  then 


^-1  =  ^n^ 


Therefore 


A"  =  An_i/cn. 


Hote:  The  first  step  of  the  above  proof  shows  that  the  matrix  AAj^_-j_  ,is  a 
scalar  matrix.  This  fact  will  be  used  later  as  a  check  for  this  method. 

Mathematical  induction  will  be  used  to  prove  that  c-^  is  eqxial  to  pj^. 
First  of  all  by  choice: 


ci   -  T(AAq)  =  T(AI)  =  T(A)  -  SI  -  pi 
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also 


2  2 

AA-j_  "  A(AAq  -c^I)    =  A     -c-j^A  =  A     -  p-^A 


so 


2c2  ■=  T(aA3_)  -  T(a  )  -  PiT(A)  -  Sg  -  P^s^  -  2p2. 


Therefore 


Cg  "  P2' 


Next 


{13  2  3  2 

A(AAq  -c-]_I)  -  C2I     "A     -  c-j_A     -  C2A  "A     -  P]_A    -  P2A 


so 


3  2 

3C3   =  T(AA2)   =  T(A  )  -  PiT(A  )   -  P2T(A)   =  s^  -  P1S2   -  P2S-L  -  3p. 


Therefore 


3"P3 


AAj^  -  A(AAj^_^  -  Cj^)   -  ArA(AAj^_2  -  Cj^_tl^)  -  Cj^l  -  A^*^  -  P3_A^  -  .   .   .  -  Pj^A 

(k+l)cj^^;j_  =  T(AAj^)   =  TCa^""^)  -  VjlU^) Pj^T(A) 

"  ^k+l  -  Pl^k  -  .    .    .  -  Pk^i  -  (k+l)Pk+i  . 


Therefore 


°k+l  "  Pk+1* 


Also  by  the  induction 


A2 


A^  =  A  -  P3_I 

2 

AAi  -  P2I  =  A    -  p^A  -  P2I 


^n  "  ^n-1  -  Pn^  "  ^"^  "  ^1^'"'^  -  .   .  .  -  P^I. 


By  the  Hamilton-Cayley  theorem  then, 


^l 
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As  an  example  of  this  method,  consider  again  the  matrix  A,  Jlrst  com- 
puting 


ci  =  1/1  T(AAq)  =  7 


h 


AA^  -  71 


-6  -1  -1  -1 

-2-3   3   0 

0  -k    -9      3 

2   10-3 


Cg  =  1/2  T(AA^)  =  -18 


12 

$ 

S 

-1 

A2  -  AA^  +  181  - 

h 
111 

23 

-13 

2U 

11 
-15 

-6 

-1 

1 

h 
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c^  =  1/3  TiAA^)   =  25 


^3  "  ^2  '  ^^"^  " 


-2$ 

-13 

-7 

-1 

31 

18 

10 

1 

-62 

-33 

-13 

-2 

U 

2 

1 

0 

AAq 


1 

0 

0 

0 

0 

-1 

0 

0 

0 

0 

-1 

0 

0 

0 

0 

-1 

c^  =  lA  T(AA3)   "  -1 


If  AA^  had  not  been  a  scalar  matrix  then  this  woiild  have  indicated  an  error 
in  the  calciilation. 
Finally, 


,-1 


A    _/c     =  A^/c, 
n-r    n        3    U 


25 

13 

7 

-3U 

-18 

-10 

62 

33 

18 

J. 

-2 

-1 

1 

2 

0 


Substituting  the  above  values  of  c,,  c«,  c-  and  c,    into  the  characteristic 
equation  gives 


|rl  -  Al   =  r' 


^  -  7r3  . 


iSr"^  -  25r  +  1  «  0. 


This  value  of  the  inverse  of  A  checks  with  results  calculated  previously. 


22 


II'JVERSION  OF  A  MTRIX  BY  PARTITIONING 
Partitioning  of  a  Matrix  to  Compute  the  Inverse 

The  motivation  for  the  method  of  inversion  by  partitioning  is  the  fact 
that  smaller  matrices  are  easier  to  invert  than  large  ones.  With  this  method, 
(2),  it  is  possible  to  invert  a  matrix  by  computing  the  inverses  of  two 
smaller  matrices  and  performing  some  matrix  multiplication  and  addition.  Two 
different  sets  of  equations  vd.ll  be  developed  which  -id-ll  each  give  the  in- 
verse .  ... 

Let  D  be  the  inverse  of  the  matrix  A,  Partition  A  between  a  pair  of 
rows  and  corresponding  columns.  Also  partition  D  and  I  similarly.  The  product 
HD  can  then  be  written  as 


■^1  ''^2 


^11,^12 
^21  °22 


I  Z 
Z  I 


Performing  the  matrix  multiplication  results  in  the  follox-jing  equality: 


%°11  ■"   ^12^21  hl^   *  ^12°22 
^1^11  ^  ^2^1   ^1^  ^  ^2^22 


I  Z 
Z  I 


In  order  for  these  two  matrices  to  be  identical,  the  corresponding  sub- 
matrices  must  be  equal;  that  is,  the  follox-ang  matrix  equations  must  hold; 


1)  ^22Pll  *  ''12^21     "^ 


2)  A-j2D2_2  ■*■  -^12^22  "   ^ 
h}     A22P22  *  -^2^22  =  ^ 
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'mJhen  the  product  AD  is  reversed  a  different  set  of  four  matrix  equations  will 
result.     The  equations  for  the  reversed  product,  DA,  are  given  below. 

^')     ^ll\l  ^  ^2^1  =  ^  2')     D^i^   .  D^i^2   =  2 

3')     D2iA^  -^  D22A21  =  Z  U')     D^^A^   +  D22A22   =  I 

To  find  D  it  is  necessary  to  solve  the  ttro  sets  of  matrix  eqxiations  for 

^11*  ■'^12*  ^21  ^"^  ^22*     ^l"^^g  ^o^  JD]_2  in  equation  2,  and  for  D21  in 
equation  3     gives 

5)     D12     =     -A1XA12D22  6)     Dg-L    -     -^22^1^11- 

Substituting  the  value  of  D-j^  in  equation  ^  into  equation  1;  one  obtains 

^l(-%^12JD22)  *  ^2^22  =  I- 
Collecting  coefficients  of  D22   gives 

^-^1%%*  ^2)^22  "  ^• 
By  definition  of  the  inverse  of  matrix 

D22  -  (A22  -  ^2.hlh2^ 


or 


7)    D22  =     (A22  -  ^j^^l^)'^. 


Solving  for  D^  in  equation  1  gives 


8)     Di3_  -     A^i  -  a{iAi2D2i. 
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Using  equations  $,   6,  7  and  8  one  can  calculate  D.  In  this  case  it  is 
necessary  to  fiixi  the  inverses  of  two  matrices,  A-^^^  and  the  inverse  of  the 
matrix  in  equation  7.  If  either  of  these  txTO  matrices  are  singular  then 
this  set  of  equations  can  not  be  used.  However,  by  using  a  different  combina- 
tion of  the  original  two  sets  of  equations,  another  set  of  foxir  equations  may 
be  obtained  in  which  two  different  inverses  must  be  found.  These  equations 
are  the  folloTdlng: 

7- )  D^j  =  Ajl  -  ql^^I,^^  B')     D^  -  (^  -  A^4^i)-^ 

In  this  set  the  inverses  of  A   and  the  matrix  in  equation  8  must  be  cal- 
culated. This  shows  that  in  order  to  use  partitioning  to  compute  the  inverse 
of  a  matrix  at  least  one  of  the  diagonal  matrices  must  be  non-singular. 

It  is  interesting  to  note  that  by  using  equations  5,  6  ,  7  and  8  it  is 
also  possible  to  calculate  D,  However,  in  this  case  the  inverses  of  all  fovir 
of  the  matrices  discussed  above  would  have  to  be  calculated.  Besides  the 
difficulty  of  calculating  two  additional  matrices,  this  set  of  equations  is 
more  restrictive  since  it  requires  two  more  matrices  to  be  non-singular. 

An  example  of  a  matrix  for  which  partitioning  would  be  advantageous  is 
a  U  X  li  matrix.  This  matrix  could  be  divided  into  fovoc   submatrices  with 
2x2  matrices  on  the  diagonal.  As  was  shoira  in  the '.second  method  of  this 
report,  the  inverse  of  a  2  x  2  matrix  is  easy  to  compute  using  the  adjoint 
matrix. 

Below  is  an  example  of  this  procedure  using  the  first  set  of  equations. 
Suppose  that  the  matrix  A  is  partitioned  so  that 


2$ 


hi   h2 
hi   h.2 


1 

-1 

1  -1 

-2 

h 

1    3 

0 

-h 

i  -2 

2 

1 

1  ° 

0 
3 
U 


-1 


Next  calculate  A-q  \ising  the  adjoint  matrix 


-1 


All  ■  VIA; 


111 


k     1 
2  1 


1/2 


i;  1 
2  1 


-1 


and  check  by  substituting  into  Ai-]_Ai]_  ■  I: 


1 

-1 

\ 

1 

1 

•1 

0 

1/2 

S3 

-2 

h 

2 

1 

0 

1 

(It  is  important  to  check  every  inverse  matrix  that  is  computed.  However,  to 
save  space  in  this  paper  the  check  will  not  always  be  shown.) 
Now  find  ' 


-1 

\     -1" 

-1      -1 

"-1  -u' 

A11A12   -  1/2 

2       1 

L              J 

3      0 

=  1/2 

1     -2 

A2lAn  "  1/2 


0    -h 

\    1" 

-1/2 

"-8 

-k 

2     1 

2    1 

10 

3 

26 


^  Ai^i2 '  ^/2 


"o    -u" 

"-1     -l" 

-2       U" 

-2       1 

1     -2 

-1/2     -5 

A22  -  ■'^i^n^ 


-2     3 

+ 

2     -li^ 

0     -1 

0   U 

X/2      5 

1/2       9 

Finally 


-1 


D22   »   (A22-  A2iAi3_Ai2)        -  2 


H                                               M 

9    1 

18     2 

1/2     0 

-1    0 

w                     « 

Next  f  iixi 


-1 


D^  =  -^^\2  "  -1/2 


-1 

'k 

18     2 

a 

7 

1 

1 

-2 

-1    0 

-10 

-1 

•   ' 

. 

. 

- 

. 

-1 


'21  =  -^22^1%  "  -^/2 


18    2 

-8 

-li" 

62    33 

-1    0 

10 

3 

-i;    -2 

^%2^21  -  V2 


-1 

-li" 

'62 

33 

=  1/2 

'-k6 

-25' 

1 

-2 

-1; 

-2 

70 

37 

^11  ■  ^  -  ^^%1-  V2 


"ll  1" 

2     1 

+  1/2 

"ii6 
-70 

25" 
-37 

s: 

"25 

-31; 

13 

-18 

27 


Therefore 


A-^ 

=  D  = 

211 
^1 

D-L2 
D22_ 

= 

25  13 

-3U  -18 

62  33 

-h  -2 


7   1 
-10  -1 


18   2 
-1   0 


J  » 


This  value  for  A~  checks  with  results  computed  previously. 

The  example  shows  how  this  method  can  be  used  to  aid  in  hand  calcula- 
tion of  the  inverses  of  small  matrices.  It  also  has  advantages  for  calcu- 
lating inverses  with  digital  computers.  Assume  that  it  is  necessary  to 
find  the  inverse  of  a  100  x  100  matrix.  The  storage  capacity  of  the  avail- 
able computer  limits  the  size  of  the  matrix  to  be  inverted  to  60  x  60  and 
the  size  of  matrices  to  be  multiplied  to  50  x  50.  By  partitioning  the  100 
X  100  matrix  into  four  50  x  50  submatrices  it  would  be  possible  to  calculate 
the  inverse  part  by  part  and  assemble  the  final  result  outside  the  computer. 
This  shows  that  by  the  use  of  partitioning  a  computer's  capacity  for  finding 
the  inverse  is  almost  doubled. 

This  method  can  be  extended  by  repeated  partitioning  to  compute  the  in- 
verses of  any  size  of  matrix.  As  an  illustration  consider  the  following 
6x6  matrix: 


-^11  "'^12 
i2i  A22 


^11  ^12  ^13  ^11; 

^21  ^22  ^23  ^21; 

"31  "32  ^33  ^3h 

%1  %2  %3  ^hh 


'51  "52  ^53  ^Sh 
^61  ^62  ^63  ^61; 


^15 

^16 

"25 

^26 

^3B 

"36 

V 

%6 

'S5 

"56 

^6^ 

^66^ 
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After  A  has  been  subdivided  in  the  above  rianner.  A,,  can  be  subdivided  into 
four  2x2  matrices.  The  inverse  of  A^  can  then  be  calculated  using  the 
method  of  partitioning.  Once  aT-,  is  known  the  inverse  of  A  can  be  computed 
using  equations  5,  6,  7  and  8.  By  using  this  process  the  inverses  of  all  the 
matrices  to  be  inverted  are  easy  to  compute.  However,  the  increase  of  the 
number  of  matrices  to  be  inverted  is  2n  where  n  is  the  number  of  times  the 
original  matrix  is  partitioned. 

The  Inversion  of  a  Matrix  by  Bordering 

As  will  be  seen  shortly,  bordering  is  very  closely  related  to  the  method 
of  partitioning  described  in  the  first  part  of  this  section.  In  this  method 
the  given  matrix  Aj,  will  be  regarded  as  the  result  of  bordering  a  matrix  of 
order  n-1.  This  submatrix  will  be  designated  by  A  .,  and  its  inverse  A~_-, 
will  be  assumed  to  be  known.  Thus 


^1 


Vl  \2n 
^In  ^nn 


^1    a^  •••  a;j_^^_^   I  a.^ 
^1    ^22  '••  ^,n-l   I  ^n 


where  ^j2n  ^^  ^  colvimn  matrix  and  Ap-j^  is  a  row  matrix. 

The  inverse  of  A^^^  will  be  designated  by  D^^  and  will  be  partitioned  in 
the  form 


D 


n 


Dn-1  Cl2n 
D21n  dn 
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where  D^-i  is  a  matrix  of  order  n-1,  D]_2n  is  a  column  matrix,  E^^^  is  a  row 
matrix,  and  d^^  is  a  number.     The  product  A^^  can  then  be  varitten  as 


Vl     %.2n 
''^In     %in 


^n-1    %n 


I     Z 
Z     I 


aJ      • 


Substituting  the  elements  into  equations  $,   6,  7  and  8  from  the  method  of 
partitioning  gives 


■"n 


Vl  *   (Vl^l2n^lnVl)^n  "Vl^ian^n 


-^InCl'^ 


dn 


where  d^^  =  (a^  -  ^iViAig)'  . 

The  inverse  of  a  matrix  A  can  be  obtained  by  using  the  above  equation  to 
compute  successively  Dj^,  D2,  D^,   .   ,   .  ,  D    where  -■.... 


-1 

%  -  [-11] 


-1 


D2 


^11  1^15 
^21  I  ^2 


^11    3-12  '  ^13 
a2i    a22  |  a23 

^31    ^32  I  ^33 


-1 
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-1 


D 


n 


^11 

a.j2 

• 

• 

• 

^l,n-l 

l-m 

^1 

^2 

• 

• 

• 

^,n-l 

^n 

^-1,1 

^-1,2 

• 

• 

• 

n 

^ 

S2 

• 

• 

• 

'^n,n-l 

Sn 

As  an  example  of  this  method  consider  the  following  matrix: 


K 


1 

-1 

-1 

-1 

2 

h 

3 

0 

0 

-h 

-2 

3 

2 

1 

0 

h 

To  compute  Dj^  using  this  method  it  would  first  be  necessary  to  find  the 
inverses  of  the  following  three  matrices: 


^1 


[^] 


A2 


•1 

1 

-1 

-1 

1 

-1 

A3- 

-2 

ii 

3 

-2 

h 

-1 

0 

-ii 

-2 

However,  since  it  is  easy  to  find  the  inverse  of  A2  using  the  adjoint  matrix 
this  will  be  done,  and  will  leave  only  two  borderings  necessary  to  find  Dr. 
Bordering  Ao  gives 
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^2    ^123 
■'^13  ^33 


1 

-1  [  -l" 
1 

2 

Ul    3 

0 

7^.^_ 

where 


Aj^  =  1/2 


h    1 
2  1 


J  • 


First  find  -aT  A223  and  ~^i^^   « 


"■^  ■^123  ~  ~  ' 


k    1 
2  1 


-1 
3 


-1/2 


-1 
1 


-^13^"''  =  -^/2[o  -i^]       =  -l/2[-8  'h] 


Next 


,  as  a  check  on  these  values,  calculate  •Aon-sC"-^  ■^123^  ^^  ^~'^13'^"'^'^123 


^13^-4^23)  =  -l/2rO    -1;]  =  2         (-^i3^'^)^23  '^  -Vsf-S  . -Ul 


-1 
3 


-1 


Since  these  valiies  are  equal,   calculate  d,     where 


S     °  ^33  "^  ^"-'^W  ■'^12^   =-2+2=0. 


Notice  from  the  equation  for  computing  D,  that  every  term  of  D^  must  be 
multiplied  by  d,,  which  in  this  case  is  undefined.  Since  this  is  impossible. 


^re^-: 


32 


the  inverse  of  Ai  cannot  be  computed  in  its  present  form  by  this  method. 
Remembering  that  two  rows  of  A  may  be  interchanged  when  finding  the  inverse 
if  after  the  inverse  is  calculated  the  corresponding  columns  of  the  inverse 
are  changed  back;  Ay   will  be  rewritten  by  interchanging  row  3  and  row  k» 


1   -1 

-1 

-1 

2       k 

3 

0 

2       1 

0 

h 

0     -ii 

-2 

3 

Since  ^2  remains  the  same,  the  first  step  is  to  border  A.^.  This  gives 


H' 


h 


^123 


^•213  ^33 


1  -1-1 
-2   lil  3 


1-- 

2   1  I  0 


As  A^     has  already  been  calculated,  solve  again  for  -A^  k-^^  and  -■^213-'^  • 


^  A3_23 


-1/2 


\  1' 

-1 

=  -1/2 

-1 

2     1 

3 

1 

-^13^^  =  -l/2[2     1 


h     1 
2     1 


-i/2rio  3.1 


Next  calculate  A2i3(-iJ  A-i_23)  and  i-i^i^k^  ^hly 
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w-'i\2i'>  -  -'4' '] 


-1 

1 


1/2 


(-^13^^  ^^3  "  -^/^[lO    3] 


-1 
3 


1/2 


-1 


Since  these  are  eqiial,  solve  for  d^    where 


d^^  -  a^^  +  (-A2i3A^'^\23^  =  0  +  1/2  =  1/2. 


Hence 


dj-a 


Now 


D2  -  Aj-  +  (^"^^3^03^^^^  "  ^/2 


U    1 
2     1 


+2 (-1/2) (-1/2) 


-1 
1 


[X0     3] 


"-6 

-2 

-3 

-1 

1/2 

m 

12 

h 

6 

2 

-1. 


^123  "  -^  ^123S  = 


1 
-1 


^213  "  -^2]3^S"  [-^°     -^] 


Therefore 


3U 


^"'  "  ^3  = 


D2        D123 
D2-L3    d3 


-3    -111 

I 
6      2   1-1 

-10     -3I  2 


Next,  repeat  this  process  after  bordering  \, 


\ 


A3  ^121; 


1  -1  -11  -1 

2  U  3  I  0 
2  1  ■  0  I| 
0  -1;  -2  j    3 


Solve  for  -aT  jUpi    and  -^2-1  kAo   . 


-A3     A;i_2i^    -    -1 


-3 

-1 

1 

-1 

-7 

6 

2 

-1 

0 

ss 

10 

10 

-3 

2 

U 

-18 

-hllih^  -  [ 


0     li     2 


-3 

-1 

1 

6 

2 

-1 

■10 

-3 

2 

[ua   0] 


-1  -1 

Check  by  shovang  that  A2]j^(-A3  A-L2i^)  and  (-A2i^A3  )A-L2ij  are  equal. 
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\ll}-'th20 


0     -h 


-] 


-7 

10 
-18 


=  "^  (-^1U^"'^-\21.  "  [i^    2     o] 


-1 
0 
ii 


-k 


Now  solve  for 


-1 


\ 


.-1 


-1 


V  ■  'kk  *  '-^ii.V,^u>  -  3  -  1.  =  -1 


d^.-l 


Next  calculate  D  ,  D  ,  and  D  ,  , 


-3 

-1 

1 

-7 

6 

2 

-1 

- 

10 

10 

-3 

2 

-18 

[k    2     O] 


25      13      1 

-3U    -13    -1 

62       33       2 


°12U  ■  -1 


.  -7 

7 

10 

B 

-10 

-18 

18 

^211" 


-iFii     2     ol  =  l-h     -2     o] 


Therefore 


=3         "l2U 


25 

13 

1  ] 

7 

•3it 

-18 

-1  1 

-10 

62 

33 

2  I 

18 

-li 

-2 

0  j 

-1 
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All  that  remains  in  calciilating  D,  is  to  interchange  column  three  and 


column  four  of  Di  .  Thus 


25 

13 

7 

1 

•31;  ■ 

-18 

-10 

-1 

62 

33 

18 

2 

-k 

-2 

-1 

0 

This  value  checks  with  results  computed  previously. 

One  advantage  of  this  method  is  the  fact  that  the  computations  may  be 
checked  when  calculating  each  d^  and  after  every  D  is  computed.  One  disad- 
vantage, as  illustrated  in  the  example  just  completed,  is  that  it  is  un- 
desirable to  have  a  submatrix  on  the  diagonal  which  does  not  have  an  inverse, 

ITERATIVE  METHODS  OF  I-ilTRIX  INVERSION  ' 


Direct  methods  of  matrix  inversion,  were  considered  in  the  three  pre- 
vious sections  of  this  report.  Since  it  is  possible  for  roundoff  errors  and 
errors  in  computation  to  occur,  a  method  of  iteration  for  improving  the 
accuracy  of  the  inverse  of  a  given  matrix  >ri.ll  be  presented.  Iterative 
methods  are  characterized  by  a  cyclicing  process  in  which  the  same  equations 
are  used  repeatedly  to  obtain  the  inverse.  An  initial  matrix  for  the  inverse 
must  be  given.  V/hen  the  desired  accuracy  is  obtained  the  process  is  stopped. 
Before  this  method  of  iteration  is  discussed,  the  basic  principle,  (U),  pp. 
177-182,  of  the  method  will  be  developed. 

This  principle  involves  the  raising  of  an  error  matrix  to  higher  and 
higher  powers.  In  order  for  the  approximations  for  the  inverse  to  approach 
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the  desired  accxiracy,  the  elements  of  these  higher  povrers  of  the  error  matrijc 
must  become  numerically  smaller  and  approach  zero.  It  will  be  shown  that  if 
the  sums  of  the  absolute  values  of  the  elements  in  each  raJ'  or  each  column  of 
a  given  matrix  are  less  than  1,  then  the  elements  of  the  higher  poi'Xers  of  the 
matrijc  approach  zero.  The  r-th  row  sxm  of  an  n  x  n  matrix,  rRSA,  is  defined 
as  this  sum  of  the  absolute  values  of  the  elements  in  the  r-th  row,  A  sim- 
pler statement  of  the  above  theorem  is,  "powers  of  a  matrix  whose  rov;  sums 
are  less  than  1  will  approach  Z,"  The  proof  to  be  given  will  consider  the 
case  of  the  "row  sums"  only.  The  proof  of  the  column  sums  would  follow  in  a 
similar  manner. 

Consider  a  matrix  B  ■vrii.ose  row  sums  are  less  than  1,  Then  there  will 
exist  at  least  one  row  sum  of  B  of  Largest  value  p  where  p<l.  It  will  be 
shown  that  the  row  sums  of  B^^  are  no  larger  than  p"^.  This  is  done  by  noting 
that  all  rov;  sums  of  a  product  of  two  matrices  AB  are  no  larger  than  p  times 
the  corresponding  row  sums  of  A,  It  will  follow  that  all  the  row  sums  of  b" 
approach  0  as  m  is  increased  because  p^  will  approach  0, 

This  result  will  be  proved  for  the  product  of  two  n  x  n  matrices  A  and 
B.  The  product  is  given  by 


,m 


n 


^  =  (^s^^^rs)  "  ^T.  a^k^^). 

k=l 


Considering  only  the  elements  in  the  r-th  row;     that  is,  the  elements  ob- 
tained for  a  fixed  r  and  for  s  varying  from  1  to  n. 


rRSAB 


n 


k=: 

n_     n 


Sk^kl 


n 


Sk^k2 


n 


k=l 


s"!  k"l  s=l  k=l 


ark 


^ks 
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< 


n 


^rk 


3=1 


^ks 


n 

z 

k=l 


'rk 


(kRSB) 


< 

>  p 

k=l 

Sk 

n 

-p> 

k=l 

^rk 

=  p(rRSA). 


This  result  shows  that  the  row  suns  of  AB  are  no  greater  than  the  product 
of  the  maxiravun  row  sum  of  B  and  the  corresponding  row  sum  of  A, 

If  the  matrix  A  =  B,  the  discussion  above  indicates  that  the  row  sums  of 

2 

B  are  less  than  or  equal  to  the  product  of  p  and  the  corresponding  row  sum 

2 

of  B.  This  means  for  this  case  that  the  row  sums  of  B  are  no  larger  than 

2  2  3 

p  .  Next  letting  A=B  ,  one  can  conclude  that  the  row  sums  of  B  are  no  larger 

2 

than  the  product  of  p  and  the  corresponding  row  sum  of  B  ,  Since  the  largest 

2  2  _3 

row  sura  of  B  is  not  greater  than  p  ,  the  row  sums  of  B^  are  all  less  than  or 

3  k 

equal  to  p  ,  Similarly,  one  can  conclude  by  letting  A  =  B  that  the  row  sums 

k+1 

of  B    are  no  larger  than  the  product  of  p  and  the  corresponding  row  sums  of 

k  k  k  k+1 

B  ,  Since  the  row  sums  of  B  are  no  larger  than  p  ,  the  row  sums  of  B    are 

k+1 
no  larger  than  p   .  Since  p  was  assumed  to  be  less  than  1,  the  powers  of  p 

must  approach  0.  Since  this  is  the  case,  there  must  exist  an  m  such  that  p 

is  approximately  equal  to  0.  This  means  that  for  this  m  the  row  sums  of  B 

must  be  approximately  equal  to  0.  Since  by  definition  the  row  sums  are  the 

sum  of  the  absolute  values  of  the  elements,  the  row  sum  approaching  0  as  the 

power  is  increased  implies  that  each  element  must  approach  0,  Therefore,  if 

the  row  sums  of  B  are  approximately  equal  to  0,  then  B  is  approximately 

equal  to  Z,  Thvis,  the  theorem  is  proven  for  the  row  sum.  As  was  stated 

earlier,  the  part  for  the  column  sums  can  be  proven  in  a  similar  manner  using 

column  sums  instead  of  row  suns. 
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Improving  the  Acciiracy  of  the  Inverse  of  a  Matrix 

With  the  background  knowledge  of  matrices  whose  powers  approach  the 
zero  matrix,  the  basic  iterative  procediire,  (li),  pp.  182-187,  for  improving 
the  accuracy  of  the  inverse  of  a  matrix  will  be  discussed.  The  first  ap- 
proximation of  the  inverse  of  A  will  be  designated  by  Dq.  This  approximation 
will  be  improved  upon  using  the  following  sequence: 

D^  -  (21  -  DqA)Dq   ■ 

Dg  =  (21  -  -^A.)!)^ 

\  =  (21  -  D^.lA)D^.l. 

It  is  obvious  that  the  error  matrix  for  Dq  can  be  written  as 

Eq  -  C  -  I  -  DqA. 

It  can  be  shown  that  the  error  matrices  of  subsequent  approximations  of  the 
inverse  are  given  by 

pin 

Therefore  if  the  first  approximation  is  such  that  the  row  stuns  or  column  suras 
of  the  error  matrix  C  are  all  less  than  1,  then  the  higher  powers  of  C  will 
approach  Z.  Since  the  error  matrix  decreases  geometrically  xjith  successive 
Dj^,  the  accuracy  should  improve  rapidly.  Once  the  elements  of  E  are  equal 
to  zero  X'/hen  rounded  to  the  desired  number  of  decimal  places  then  D_  will  be 


Uo 


equal  to  A~  to  the  desired  accuracy. 

As  an  example  of  this  method  consider  .the  follov;ing  matrix: 


B 


1  -2  -1 
U  -3  1 
3  -1  -2 


J  • 


As  a  first  approximation  of  the  inverse,  let 


-0.3  0.2  0.2 
-0.5  0.0  0.1 
-0.2  0.3  -0.3 


Form  the  product 


DoB  = 


-0.3  0.2  0.2 
-0.5  0.0  0.1 
-0.2  0.3  -0.3 


1  -2  -1 
i;  -3   1  = 
3  -1  -2 


1.1  -0.2  0.1 

-0.2   0.9  0.3 

0.1  -0.2  1.1 


Subtracting  this  product  from  the  identity  matrix  gives  the  error  matrix 


I-DqB 


10     0 

0     10 

- 

0     0     1 

J 

1.1  -0.2  0.1 

-0.2   0.9  ,0.3 

0,1  -0.2  1.1 


-0.1  0.2  -0.1 

0.2  0.1  -0.3 

-0.1  0.2  -0.1 


J  • 


la 


Since  the  row  sums  of  C  are  all  less  than  1,  this  value  of  Dq  will  work. 
If  six  place  accuracy  is  desired  compute  certain  powers  of  C  to  determine  how 
many  iterations  will  be  necessary.  These  powers  will  also  serve  as  a  check  on 
the  computations.  The  powers  needed  are 


0.06  -0.02  -O.Oii 
0.03  -0.01  -0.02 
0.06  -0.02  -O.OU 


c^  = 


0.0006  -0.0002  -o.oooU 
0.0003  -0.0001  -0.0002 

0.0006  -0.0002  -O.oooU 


0.00000006  -0.00000002  -o.oooooooU 
0.00000003  -0.00000001  -0.00000002 
0.00000006  -0.00000002  -O.OOOOOOOli 


J  • 


8  -1 

Since  C  equals  Z  to  six  places  Do  will  equal  B   to  the  same  accuracy.  To 

compute  Do  it  is  first  necessary  to  compute  D^^. 


D-L  -  (21  -  DqB)Dq  - 


0.9  0.2  -0.1 

0.2  1.1  -0.3 

-0.1  0.2   0.9 


-0.3  0.2  0.2 
-0.5  0.0  0.1 
-0.2  0.3  -0.3 


and 


D^B 


-0.35   0.15   0.23 

-0.55  -0.05  0.2U 

-0.25   0.25  -0.27 


1  -2  -1 
U  -3   1  = 
3  -1  -2 


0.9U  0.02  O.Oli 
-0.03  1.01  0.02 
-0.06  0.02  I.OU 


hz 


As  a  check  D,B  subtracted  from  the  identity  matrix  should  equal  C  •  It  is 
easy  to  see  that  this  is  true. 

The  computation  for  D^  and  D^B  gives 


D2  «  (21  -  D^B)D^ 


1.06  -0.02  -O.Oii 
0.03  0.99  -0.02 
0.06  -0.02   0.96 


-0.35  0.15  0.23 

-0.55  -0.05   0.21; 
-0.25   0.25  -0.27 


and 


DgB 


-0.3500  0.1500  0.2U98 
-0.5500  -0.0500  0.2U99 
-0.2500   0.2500  -0.2502 


1  -2  -1 

h    -3   1 
3  -1  -2 


0.999U  0.0002  O.OOOl; 
-0.0003  1.0001  0.0002 
-0.0006  0.0002  l.OOOU 


In  comparing  this  with  the  identity  matrix  one  sees  that  the  difference  is 
equal  to  C  .  Finally 


D^  =  (21  -  D2B)D2 


1.0006  -0.0002  -O.OOOU 
0.0003  0,9999  -0.0002 
0.0006  -0.0002   0.9996 


-0.3500  0.1500  0.2^98 
-0.5500  -0.0500  0.2lt99 
-0.2500   0.2500  -0.2502 


-0.35000000  0.15000000  0.2U999998 
-0.55000000  -0.05000000  Q.2h999999 
-0.25000000   0.25000000  -0.25000002 


U3 


If  the  product  D^B  is.  formed  and  subtracted  from  the  identity  matrix,  the 

Q 

result  will  equal  C  .  This  is  a  final  check  of  computational  errors. 
Since  six  place  accuracy  was  desired  Do  can  be  rounded  off  sis 


"3 


-0.350000  0.1^0000  0.250000 
.-0.550000  -0.050000  0.250000 
-0.250000  0.250000  -0.250000 


In  this  case  Do  is  the  exact  inverse  of  B. 

It  is  interestir-g  to  note  what  would  happen  if  some  Dj^  is  the  exact 
solution  so  that  Dj^  =  I.  If  this  is  the  case  then, 

^k+1  =  (21  -  Dj^)Dj^  =  (21  -  1)D^   =  IDk  =  Dk- 

Thus  one  criterion  for  halting  the  iteration  process  is  that  two  subsequent 
approximations  are  equal.  Another  process  for  determining  the  number  of 
iterations  needed  is  the  one  used  in  the  example.  This  was  computing  the 
powers  of  C  until  an  error  matrix  vrLth  small  enough  elements  is  obtained. 

Inversion  Method  for  Strongly  Diagonal  tiatrices 

For  certain  matrices  it  is  easy  to  guess  an  initial  approximation  of  the 
inverse  which  leaves  the  row  siims  of  the  error  matrix  C,  less  than  1,  •  These 
matrices,  (U),  pp.  187-189,  are  called  strongly  diagonal  matrices.  Tciey  are 
characterized  by  diagonal  elements  whose  absolute  values  are  larger  than  the 
sums  of  the  absolute  values  of  the  remaining  elements  of  each  of  their  rows, 
that  is,  iRSA  <  2|a^^  I  for  each  row.  If  each  of  the  rows  of  a  strongly 
diagonal  matrix  is  divided  by  the  diagonal  element  then  the  row  sums  of  the 


hh 


new  matrix  are  less  than  2,  If  this  matrix  is  subtracted  from  the  identity- 
matrix  then  the  row  sums  of  the  resulting  matrix  are  less  than  1,  Remembering 
the  equation  for  the  error  matrix  is  I  -  DA  it  follows  that  the  matrix  to 
pick  for  D  is 


^11 


0  .  .  . 


0  3^2  .  .  . 


0  . 


nn 


This  value  for  Dq  makes  C  the  matrix  described  in  the  above  paragraph  whose 
row  sums  were  less  than  one. 

As  an  example  of  a  strongly  diagonal  matrix  consider  the  following  matrix: 


B  = 


h 

-2 

-1 

3 

-5 

1 

-2 

h 

8 

First  checking  the  rows  gives 


row  1:  1  +  2  =  3  <  ii 
row  2:  3+l  =  i|<5 
row  3:     2+U  =  6<8 


which  verifies  that  matrix  B  is  indeed  a  strongly  diagonal  matrix.     Let 


h$ 


^0  = 


0.25  0  0 
0  -0.2  0 
0     0    0.125 


Then 


C  =  I  -  DqA 


1 

0 

0 

m 

0 

1 

0 

- 

0 

0 

1 

1.0   -0.5  -0.25 
-0.6    1,0  -0.2 

-0.25  0.5  1.0 


0-    0.5  0.25 

0.6    0    0.2 

0.25  -0.5  0 


Checking  the  row  sums  of  C 


IRSC  =  0.50  +  0.25  =  0.75 
2RSC  =  0.60  +  0.20  =  0.80  . 
3RSG  =  0.25  +  0.50  =  0.75  . 


Since  for  this  value  of  D-  the  error  matrix  C  satisfies  the  property  described 
in  the  introduction  of  this  section,  B  "^  can  be  found  using  the  iterative 
method  discussed  earlier  in  this  section. 

This  method  can  be  used  for  calculating  the  inverse  of  a  matrix -that  is 
not  strongly  diagonal  if  the  matrix  can  be  transformed  into  a  strongly  diagonal 
matrix  by  interchanging  rows.  The  first  step  is  to  place  the  element  of 
largest  absolute  value  on  the  diagonal.  This  can  be  done  by  interchanging 
two  rows.  Next,  eliminate  the  row  and  column  where  this  element  is  contained, 
and  place  the  element  of  largest  absolute  value  in  the  remaining  submatrix  on 
the  diagonal.  Continue  this  process  until  the  size  of  the  submatrix  is  one. 


k6 


Next,  check  the  row  sums  of  the  newly  formed  matrix  to  determine  whether  the 
matrix  is  strongly  diagonal.  If  so,  then  its  inverse  can  be  computed  using, 
the  procedure  discussed  in  this  part.  To  find  the  inverse  of  the  original 
matrix,  it  is  necessary  to  interchange  the  corresponding  columns  of  the  cal- 
culated inverse  in  reverse  order. 

Instead  of  interchanging  rows,  it  would  work  equally  well  if  the  columns 
of  the  original  matrix  were  interchanged.  Then  the  corresponding  rox^s  of  the 
calculated  inverse  XTOuld  be  interchanged  in  reverse  order.  An  obvious  disad- 
vantage of  this  procedure  is  the  necessity  of  remembering  the  interchange  and 
their  order. 

CONCLUSION 

Only  a  few  of  the  many  different  methods  of  computing  the  inverse  of  a 
matrix  are  given  in  this  report.  Some  of  the  methods  not  included  are  varia- 
tions of  those  contained  here  and  others  are  methods  which  can  only  be  used  on 
special  kinds  of  matrices. 

Solving  the  system  of  linear  equations 

AX  =  G 

using  the  inverse  of  the  coefficient  matrix  A  is  only  one  technique  of  solu- 
tion. It  becomes  advantageous  to  use  this  technique  only  when  the  coefficient 
matrix  is  large  or  when  the  same  coefficient  matrix  A  is  used  for  several  dif- 
ferent constant  matrices  G, 
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With  the  advent  of  high-speed  electronic  digital  computers  it  became 
feasible  to  soIvb  large  systenis  of  equations.  Many  techniques  have  been 
developed  to  solve  these  equations.  One  of  these  involves  calculating  the 
inverse  of  the  coefficient  matrix. 

Using  matrix  notation,  a  system  of  n  linear  equations  in  n  unknowns 
can  be  ^•nritten  in  the  matrix  form:       •     '  - 

AX  =  G. 


-1 


By  multiplying  both  sides  of  the  equation  by  A 

A~\x  =  A'h, 
X  =  A"-^. 


Thus,  if  the  value  of  the  inverse  of  the  coefficient  matrix  is  knovxn,  the 
solution  for  the  system  can  be  found  by  one  matrix  multiplication. 

The  purpose  of  this  report  is  to  consider  eight  different  methods  by 
which  the  inverse  of  the  coefficient  matrix  may  be  calculated.  For  each 
method  an  effort  is  made  to  explain  the  method,  prove  its  validity,  and  givB 
an  example  using  the  method.  Also,  for  most  of  the  methods  advantages  aM  ' 
disadvantages  are  given.  The  eight  different  methods  of  matrix  inversion 
are  discussed  in  four  sections.  The  first  section  considers  matrix  inver- 
sion using  elementary  operations.  Of  all  the  methods  discussed  in  this 
paper,  those  methods  which  involve  elementary  operations  are  the  most  easily 
adapted  to  machine  computation.  Only  the  most  general  case  of  these  methods 
is  discussed  here. 


The  second  section  considers  methods  of  matrix  inversion  involving 
the  determinant  oi  a  matrix.  The  three  methods  incliided  use  the  adjoint 
matrix,  the  characteristic  equation,  and  the  trace  to  calcvdate  the  inverse. 

The  inversion  of  a  matrix  by  partitioning  and  a  special  case  in. which 
the  niatrix  is  partitioned  by  bordoring  ara  developed  in  the  third  Gsetion. 

Iterative  methods  of  matrix  inversion  are  considered  in  the  fourth 
section.  These  methods  are  important  for  they  can  be  used  to  improve  the 
accuracy  of  inverses  calculated  by  other  methods,  A  basic  method  of  itera- 
tion and  a  particular  type  of  matrix  for  i-rhich  a  good  starting  value  is 
easily  determined  are  discussed. 

Only  a  few  of  the  many  different  methods  of  computing  the  inverse  of 
a  matrix  are  contained  in  this  report.  Some  of  the  methods  not  included 
are  variations  of  those  contained  here  ssd   others  are  methods  which  can 
only  be  used  on  special  kinds  of  matrices. 

Solving  a  system  of  linear  equations  using  the  inverse  of  the  coeffi- 
cient matrix  is  only  one  technique  of  solution.  It  becomes  advantageous 
only  when  the  coefficient  matrix  is  large  or  when  the  same  coefficient 
matrix  A  is  used  for  several  different  constant  matrices  G, 


